library(ggplot2)
library(tidyr)
library(tidyverse)
library(estimatr)
library(coefplot)
library(patchwork)

# Figure 5: Supporting Nuclear Use At A Specific Threshold + Average Treatment Effect

setwd("/Users/dmin/Dropbox/Replication_Data/Paper_3/")

Dat <- read.csv("Paper_3_Data_Clean.csv", stringsAsFactors = FALSE)

Dat <-
  Dat %>% 
  mutate(
    By30 = as.numeric(outcome_custom<=30),
    By50 = as.numeric(outcome_custom<=50),
    By60 = as.numeric(outcome_custom<=60),
    By80 = as.numeric(outcome_custom<=80),
    By90 = as.numeric(outcome_custom<=90),
    By99 = as.numeric(outcome_custom<=99)
  )

gg_df_1 <-
  Dat %>% 
  group_by(Treatment, NuclearTheory) %>% 
  summarize(tidy(lm_robust(By99 ~ 1, data = cur_data())))

gg_df_1 <-
  gg_df_1 %>% 
  mutate(term_label =
           factor(
             Treatment,
             levels = c("Control",
                        "FriendPro-NuclearMoral",
                        "FriendPro-NuclearUtilitarian",
                        "SenatorPro-NuclearMoral",
                        "SenatorPro-NuclearUtilitarian"),
             labels = c(
               "Control",
               "Peer\nMoral",
               "Peer\nUtilitarian",
               "Elite\nMoral",
               "Elite\nUtilitarian"
             )
           ),
         
         theory_label = factor(NuclearTheory, levels = c(1, 2, 3, 4),
                               labels = c("Never-Ever", "Deterrence", "Offense/Defense","Utilitarian"))
  )




label_df <- 
  gg_df_1 %>% 
  filter(term_label == "Elite\nUtilitarian")

pos <- position_dodgev(height = 0.5)

g1 <- 
  ggplot(gg_df_1, aes(x = estimate, y = term_label, color = theory_label, group = theory_label,shape= theory_label)) +
  geom_point(position = pos) +
  geom_linerange(aes(xmin = conf.low, xmax = conf.high), position = pos) +
  geom_text(data = label_df, aes(label = theory_label, x = conf.low - 0.06), position = position_dodgev(height = 0.5)) +
  theme_bw() + 
  scale_color_manual(values=c("green3","yellow3","orange3", "red3")) + 
  theme(legend.position="none") +
  ylab("Treatment")+
  xlab("Proportion Supporting Nuclear Use at or Below 99% Chance of Success")

g1



gg_df_2 <-
  Dat %>% 
  group_by(NuclearTheory) %>% 
  summarize(tidy(lm_robust(By99 ~ Treatment, data = cur_data()))) %>% 
  filter(term != "(Intercept)")

gg_df_2 <-
  gg_df_2 %>% 
  mutate(term_label =
           factor(
             term,
             levels = c(
               "TreatmentFriendPro-NuclearMoral",
               "TreatmentFriendPro-NuclearUtilitarian",
               "TreatmentSenatorPro-NuclearMoral",
               "TreatmentSenatorPro-NuclearUtilitarian"),
             labels = c(
               "Peer\nMoral",
               "Peer\nUtilitarian",
               "Elite\nMoral",
               "Elite\nUtilitarian"
             )
           ),
         
         theory_label = factor(NuclearTheory, levels = c(1, 2, 3, 4),
                               labels = c("Never-Ever", "Deterrence", "Offense/Defense","Utilitarian"))
  )


label_df_2 <- 
  gg_df_2 %>% 
  filter(term_label == "Elite\nUtilitarian")

pos <- position_dodgev(height = 0.5)

g2 <- 
  ggplot(gg_df_2, aes(x = estimate, y = term_label, color = theory_label, group = theory_label, shape= theory_label)) +
  geom_point(position = pos) +
  geom_linerange(aes(xmin = conf.low, xmax = conf.high), position = pos) +
  geom_text(data = label_df_2, aes(label = theory_label, x = conf.low - 0.05), position = position_dodgev(height = 0.5)) +
  geom_vline(xintercept=0, linetype="dotted")+
  theme_bw() + 
  scale_color_manual(values=c("green3","yellow3","orange3", "red3")) + 
  theme(legend.position="none") +
  ylab("")+
  xlab("Average Treatment Effect Estimate")+
  coord_cartesian(ylim=c(0,4))


g <- g1 + g2

g

ggsave("Image_5.png", plot=g, width=4000, height=2200, units=c("px"))


